Metodo delle potenze inverse

Nell'analisi numerica, il metodo delle potenze inverse è un algoritmo iterativo per il calcolo degli autovettori di una matrice. L'algoritmo permette di stimare un autovettore quando è già conosciuta una approssimazione dell'autovalore corrispondente. Questo metodo è concettualmente simile al metodo delle potenze e nacque per calcolare le frequenze di risonanza nel campo della meccanica strutturale. ^[1]

Il metodo delle potenze inverse inizia con una approssimazione $\mu$ per l'autovalore corrispondente all'autovettore desiderato e un vettore $b_{0}$ , scelto sia casualmente oppure da un'approssimazione dell'autovettore. L'algoritmo è descritto dall'iterazione

b_{k+1}={\frac {(A-\mu I)^{-1}b_{k}}{C_{k}}},

dove $C_{k}$ sono qualche costanti di solito scelte come $C_{k}=\|(A-\mu I)^{-1}b_{k}\|$ . Poiché gli autovettori sono definiti a meno di uno scalare moltiplicativo, la scelta di $C_{k}$ può essere arbitraria nella teoria; aspetti pratici della scelta di $C_{k}$ sono discussi sotto.

Ad ogni iterazione, il vettore $b_{k}$ è moltiplicato per l'inversa della matrice $(A-\mu I)^{-1}$ e normalizzata. È esattamente la stessa formula del metodo delle potenze, eccetto la sostituzione di $A$ con $(A-\mu I)^{-1}$ . Più l'approssimazione $\mu$ è scelta vicino all'autovalore, più velocemente l'algoritmo converge; tuttavia, scelte sbagliate di $\mu$ possono portare a una convergenza lenta oppure a un altro autovettore rispetto a quello desiderato. Nella pratica, questo metodo è utilizzato quando si conosce una buona approssimazione dell'autovalore, e quindi servono solo poche (molto spesso anche solo una) iterazioni.

Teoria e convergenza del metodo

L'idea base del metodo delle potenze è scegliere un vettore iniziale $b$ (l'autovettore approssimato o un vettore casuale) e iterativamente calcolare $Ab,A^{2}b,A^{3}b,...$ . eccetto per un insieme di misura nulla, per ogni vettore iniziale, il risultato convergerà a un autovettore corrispondente all'autovalore dominante.

Il metodo delle potenze inverse fa lo stesso con la matrice $(A-\mu I)^{-1}$ , perciò converge all'autovettore corrispondente all'autovalore dominante della matrice $(A-\mu I)^{-1}$ . Gli autovalori di questa matrice sono $(\lambda _{1}-\mu )^{-1}$ , $\ldots$ , $(\lambda _{n}-\mu )^{-1}$ , dove $\lambda _{i}$ sono gli autovalori di $A$ . Il più grande di questi valori corrisponde al più piccolo di $(\lambda _{1}-\mu )$ , $\ldots$ , $(\lambda _{n}-\mu )$ . Gli autovettori di $A$ e di $(A-\mu I)^{-1}$ sono i soliti, poiché

Av=\lambda v\,\Leftrightarrow \,(A-\mu I)v=\lambda v-\mu v\,\Leftrightarrow \,(\lambda -\mu )^{-1}v=(A-\mu I)^{-1}v

Conclusione: Il metodo converge all'autovettore della matrice $A$ corrispondente all'autovalore più vicino a $\mu$ .

In particolare, prendendo $\mu =0$ si osserva che $(A)^{-k}b$ converge all'autovettore corrispondente all'autovalore di $A$ con il valore assoluto più piccolo.

Velocità di convergenza

Si analizzi la velocità di convergenza del metodo.

Il metodo delle potenze converge linearmente al limite, più precisamente:

$\mathrm {Distanza} (b^{\mathrm {ideale} },b_{\mathrm {Metodo~Potenze} }^{k})=O\left(\left|{\frac {\lambda _{\mathrm {sottodominante} }}{\lambda _{\mathrm {dominante} }}}\right|^{k}\right),$

quindi per il metodo delle potenze inverse si hanno risultati simili:

$\mathrm {Distanza} (b^{\mathrm {ideale} },b_{\mathrm {Potenze~Inverse} }^{k})=O\left(\left|{\frac {\mu -\lambda _{\mathrm {piu'~vicino~a~} \mu }}{\mu -\lambda _{\mathrm {secondo~piu'~vicino~a~} \mu }}}\right|^{k}\right).$

Questa è la formula chiave per comprendere la convergenza dell'algoritmo. Si osserva che se $\mu$ è scelto abbastanza vicino a un qualche autovalore $\lambda$ , per esempio $\mu -\lambda =\epsilon$ , ogni iterazione migliorerà la precisione di $|\epsilon |/|\lambda +\epsilon -\lambda _{\mathrm {piu'~vicino~a~} \lambda }|$ volte. (Si è usato che per $\epsilon$ abbastanza piccoli, "più vicino a $\mu$ " e "più vicino a $\lambda$ " sono la stessa cosa.) Per $|\epsilon |$ abbastanza piccoli, è approssimativamente uguale a $|\epsilon |/|\lambda -\lambda _{\mathrm {piu'~vicino~a~} \lambda }|$ . Dunque se si è in grado di trovare $\mu$ tale che $\epsilon$ sarà abbastanza piccolo, allora saranno sufficienti poche interazioni.

Complessità

L'algoritmo richiede di risolvere un sistema lineare o di calcolare l'inversa di una matrice. Per matrici non strutturate (non sparse, di Toeplitz,...) questo richiede $O(n^{3})$ operazioni.

Opzioni per l'implementazione

Il metodo è definito dalla formula:

b_{k+1}={\frac {(A-\mu I)^{-1}b_{k}}{C_{k}}},

Ci sono, tuttavia, multiple opzioni per la sua implementazione.

Calcolare l'inversa della matrice o risolvere un sistema lineare

Si può riscrivere la formula nella seguente maniera:

(A-\mu I)b_{k+1}={\frac {b_{k}}{C_{k}}},

sottolineando che per trovare la successiva approssimazione $b_{k+1}$ si può risolvere un sistema di equazioni lineari. Ci sono due opzioni: si può scegliere un algoritmo che risolve un sistema lineare, oppure si può calcolare l'inversa $(A-\mu I)^{-1}$ e dopo applicarla al vettore. Entrambe le soluzioni hanno complessità $O(n^{3})$ , il numero esatto dipende dal metodo scelto.

La scelta dipende anche dal numero di iterazioni. Ingenuamente, se a ogni iterazione si risolve un sistema lineare, la complessità sarà $K\cdot O(n^{3})$ , dove $K$ è il numero di iterazioni; similmente, calcolare la matrice inversa e applicarla ad ogni iterazione ha complessità $K\cdot O(n^{3})$ . Si osserva, tuttavia, che se la stima dell'autovalore $\mu$ rimane costante, per entrambi i metodi si può ridurre la complessità a $O(n^{3})+K\cdot O(n^{2})$ . Calcolare l'inversa una volta memorizzandola per applicarla ad ogni iterazione è di complessità $O(n^{3})+K\cdot O(n^{2})$ . Immagazzinare una decomposizione di $(A-\mu I)$ e usare il metodo di sostituzione per risolvere il sistema di equazioni è ancora di complessità $O(n^{3})+K\cdot O(n^{2})$ .

Invertire la matrice ha tipicamente un maggiore costo iniziale, ma minore ad ogni iterazione. Al contrario, risolvere il sistema di equazioni lineari ha un costo iniziale minore, ma richiede più operazioni a ogni iterazione.

Tridiagonalizzazione, Forma di Hessenberg

Se è necessario effettuare molte iterazioni (opoche, ma per molti autovettori), allora può essere saggio portare la matrice prima nella forma di Hessenberg superiore (per matrice simmetriche è la forma tridiagonale), che costa ${\begin{matrix}{\frac {10}{3}}\end{matrix}}n^{3}+O(n^{2})$ operazioni aritmetiche usando una tecnica basata sulla trasformazione di Householder, con una sequenza finita di similitudini ortogonali, come fosse una decomposizione QR da entrambi i lati.^[2]^[3] (Per decomposizioni QR, le rotazioni di Householder sono moltiplicate solo a sinistra, ma nel caso della forma di Hessenberg si moltiplica sia a destra che a sinistra). Usando sempre questa tecnica, per le matrici simmetriche questa procedura costa ${\begin{matrix}{\frac {4}{3}}\end{matrix}}n^{3}+O(n^{2})$ operazioni.^[2]^[3]

Per le matrici tridiagonali, La soluzione al sistema lineare per le matrici costa $O(n)$ operazioni, quindi la complessità cresce come $O(n^{3})+K\cdot O(n)$ , dove $K$ è il numero di applicazioni, che è migliore dell'inversione diretta. Tuttavia, per poche iterazioni una tale trasformazione è poco pratica.

Inoltre trasformazioni nella forma di Hessenberg coinvolge radici quadrate e divisioni, che non sono universalmente supportate dall'hardware.

Scelta della costante di normalizzazioneC_k

Nei processori General purpose (per esempio prodotti da Intel), il tempo di esecuzione di addizione, moltiplicazione e divisione è circa lo stesso. Tuttavia su hardware integrati e/o a basso consumo (DSP, FPGA, ASIC), la divisione può non essere supportata così dovrebbe essere evitata. La scelta $C_{k}=2^{n_{k}}$ ,permette divisioni veloci senza supporti hardware espliciti, dal momento che la divisione per una potenza di due può essere implementata sia come uno spostamento bit a bit (per l'aritmetica in virgola fissa) sia con sottrazione di $k$ dall'esponente (per l'aritmetica in virgola mobile).

Quando viene implementato l'algoritmo con la rappresentazione in virgola fissa, la scelta della costate $C_{k}$ è particolarmente importante. Piccoli valori porteranno ad una crescita veloce della norma di $b_{k}$ e quindi all'overflow; grandi valori di $C_{k}$ causerà l'avvicinarsi di $b_{k}$ a zero.

Applicazioni

La principale applicazione del metodo è la situazione in cui si possiede un'approssimazione di un autovalore e si deve trovare la stima del corrispondente autovettore. In tale situazione il metodo delle potenze inverse è il principale algoritmo da usare.

Metodi per trovare autovalori approssimati

Tipicamente questo metodo è usato in combinazione con altri per approssimare autovalori: l'esempio standard è l'algoritmo di bisezione per gli autovalori; un altro esempio è il metodo del quoziete di Rayleigh, che è in realtà è il metodo delle potenze inverse con l'utilizzo del quoziente di Rayleigh del vettore precedente come approssimazione dell'autovalore.

Ci sono alcune situazioni in cui si possono usare solo determinati metodi, tuttavia sono alquanto marginali.

Norma della matrice come approssimazione dell'autovalore dominante

Si può stimare facilmente l'autovalore dominante per ogni matrice. Per ogni norma indotta, è vero che $\left\|A\right\|\geq |\lambda |,$ per ogni autovalore $\lambda$ . Perciò prendendo la norma della matrice come approssimazione dell'autovalore, si può vedere che il metodo converge all'autovettore dominante.

Stime statistiche

In alcune applicazioni in tempo reale, si ha bisogno di trovare autovettori con una velocità di milioni di matrici al secondo. In tali applicazioni, tipicamente la statistica delle matrici è conosciuta a priori e si può prendere come autovalore approssimato il valore medio di quest'ultimo per qualche grande matrice campione. Ancora meglio, si può calcolare la media dei rapporti tra gli autovalori e la traccia o la norma della matrice e stimare l'autovalore medio grazie al valore medio di questo rapporto. Chiaramente questo metodo può essere usato solo con discrezione e solo quando non è cruciale una altissima precisione. Questo approccio di determinare un autovalore medio può essere combinato con uno degli altri metodi per evitare errori eccessivamente grandi.

Note

^ Ernst Pohlhausen, Berechnung der Eigenschwingungen statisch-bestimmter Fachwerke, ZAMM - Zeitschrift für Angewandte Mathematik und Mechanik 1, 28-42 (1921).
^ ^a ^b James W. Demmel, Applied Numerical Linear Algebra, Philadelphia, PA, Society for Industrial and Applied Mathematics, 1997, ISBN 0-89871-389-7, MR 1463942.
^ ^a ^b Lloyd N. Trefethen and David Bau, Numerical Linear Algebra (SIAM, 1997).

Voci correlate

Collegamenti esterni

Metodo delle potenze inverse per trovare autovettori, physics.arizona.edu

Portale Matematica: accedi alle voci di Wikipedia che trattano di matematica

[Pohlhausen-1] Ernst Pohlhausen, Berechnung der Eigenschwingungen statisch-bestimmter Fachwerke, ZAMM - Zeitschrift für Angewandte Mathematik und Mechanik 1, 28-42 (1921).

[Demmel-2] James W. Demmel, Applied Numerical Linear Algebra, Philadelphia, PA, Society for Industrial and Applied Mathematics, 1997, ISBN 0-89871-389-7, MR 1463942.

[Trefethen-3] Lloyd N. Trefethen and David Bau, Numerical Linear Algebra (SIAM, 1997).

[1]

[2]

[3]

Metodo delle potenze inverse

Indice